Journal of Chemical Theory and Computation
● American Chemical Society (ACS)
Preprints posted in the last 90 days, ranked by how well they match Journal of Chemical Theory and Computation's content profile, based on 140 papers previously published here. The average preprint has a 0.08% match score for this journal, so anything above that is already an above-average fit.
Elgendy, A.; Zeipelt, A. P.; Schäfer, L. V.
Show abstract
Molecular dynamics (MD) simulations of slow biomolecular processes, such as exploration of the conformational ensembles of intrinsically disordered proteins (IDPs), are computationally demanding. Although coarse-grained (CG) models can substantially speed up the simulations compared to all-atom MD, the sampling challenge can still be significant for large systems and long time scales. Here, we present light Martini water, a low-viscosity water model that accelerates sampling in MD simulations with the Martini CG force field. We systematically reduced the mass of the Martini water beads and verified stable, accurate integration of the equations of motion with 20 fs time steps, as typically used in Martini simulations. Light Martini water has a reduced mass of 20 amu (compared to 72 amu in the standard water model), yielding up to a 2.68-fold increase in the sampling rate of IDP chain reconfiguration in water and a 16 % increase in the lateral diffusion of lipids in a POPC bilayer. Equilibrium properties remained unaffected by the mass scaling, and the speedup was achieved without compromising simulation accuracy. The water model is trivial to implement, has no computational overhead, and should be universally applicable to Martini simulations.
Kozlowski, N.; Nagy, G.; Ottoway, C. F.; Bock, L. V.; Vaiana, S. M.; Grubmuller, H.
Show abstract
Atomistic simulations of intrinsically disordered proteins (IDPs) are notoriously sensitive to force field inaccuracies, either regarding the protein or the water model, yielding inaccurate observables such as compactness, secondary structure propensities, or kinetics. The currently most widely used IDP force fields are Amber99sb-disp (A99disp) and Charmm36m (C36m). A99disp includes a new water model and thus optimized both, the protein and the water interactions. In contrast, C36m used the Tip3p water model and optimized only protein interactions. In many cases, C36m+Tip3p underestimates radii of gyration compared to FRET or SAXS experiments. Such overly compact structural ensembles are believed to arise from an imbalance between protein-protein, protein-water, and water-water interactions, which might be due to Tip3p inaccuracies. Here, we aim at re-balancing these interactions by combining C36m with the Optimal Point Charge (OPC) water model. Recently, this C36m+OPC combination showed improved accuracy for the disordered domain of the measles virus nucleoprotein. Here we present a systematic assessment, comparing C36m+OPC, C36m+Tip3p, C36m+Tip4p, C22*, A03ws, A99sb-ws, and A99disp for five IDPs, as well as a subset of those for five globular proteins, a set of disordered AGQ-repeat peptides, and the fast folding miniprotein CLN025. We compared extensive MD simulations (> 8.5 ms) with SAXS, NMR, circular dichroism, photo-induced electron transfer (PET), T-jump infrared spectroscopy, and X-ray crystallography measurements. We found that combining C36m with OPC improved accuracy over C36m+Tip3p for IDP ensembles without compromising its accuracy for globular proteins. While also the kinetics of the AGQ-peptides were more accurate for C36m+OPC, those of CLN025 folding were less accurate. Overall, C36m+OPC showed similar accuracy as A99sb-ws and A99disp, the latter is currently considered among the most accurate protein force fields.
Akinyemi, O.; Kierzek, E.; Kierzek, R.; McSally, J. P.; Puthenpeedikakkal, A. M. K.; Mathews, D. H.
Show abstract
Umbrella sampling is widely used to estimate biomolecular free energy landscapes and relative folding stabilities. Although equilibration is a critical component of umbrella sampling workflows, the impact of equilibration length on thermodynamic predictions remains poorly understood. Here, we investigate the effect of equilibration length on relative folding free energy predictions for four RNA hairpins with loop sequences GUGAAA, CUGGGA, GUAAUA, and UUAAUU with helical stems of three base pairs. Umbrella sampling simulations were performed using an end-to-end distance reaction coordinate spanning 15-45 [A], where equilibrium simulations (windows) were spaced at roughly 1 [A] intervals. In these calculations, the hairpin stem-loops were allowed to equilibrate in an end-to-end distance window and then the coordinates were transferred to the next larger end-to-end distance window to equilibrate. Two equilibration lengths, 2 ns and 100 ns per window, were followed by 600 ns of production sampling. Potential of mean force (PMF) profiles were reconstructed using the Weighted Histogram Analysis Method (WHAM) and used to calculate pairwise free energy differences with thermodynamic cycles. Increasing the equilibration length produced substantial, sequence-dependent changes in the reconstructed free energy landscapes. The 100 ns protocol generated markedly flatter PMFs for GUGAAA and UUAAUU and pronounced reshaping of the free energy landscape for GUAAUA. These changes were accompanied by reductions in hydrogen-bonding and stacking interactions, particularly within the intermediate regions of the reaction coordinate. The resulting thermodynamic predictions, as free energy change differences, were therefore highly sensitive to equilibration length. Across nearly all hairpin pairs, the 100 ns equilibration yielded substantially larger magnitude free energy change difference values than the corresponding 2 ns equilibration, with differences that greatly exceeded replica-to-replica variability. Comparison with optical melting measurements and nearest-neighbor thermodynamic predictions revealed that 2 ns equilibration times more closely agreed with experimental values than those obtained using 100 ns equilibration. These findings demonstrate that longer equilibration can systematically alter the structural ensembles sampled during umbrella sampling and amplify predicted stability differences without improving agreement with experiment. More widely, our results highlighted equilibration length as a critical and nontrivial parameter in RNA free energy calculations and demonstrate that increased equilibration does not necessarily lead to more accurate thermodynamic predictions.
Scaramozzino, D.; Ruggeri, E.; Feshyna, M.; Lee, B. H.; Orellana, L.
Show abstract
Protein dynamics is fundamental to understand mechanisms. Although atomistic Molecular Dynamics (MD) remains the gold standard for predicting protein motions, its computational cost limits applications at large scales. Here, we present a coarse-grained (CG) simulation framework that combines Elastic Network Models (ENMs) with implicit-solvent Brownian Dynamics (BD) to efficiently generate protein conformational ensembles from native states. Benchmarking the method on more than 8,500 proteins against large datasets of atomistic MD trajectories and multi-state ensembles from Nucleic Magnetic Resonance (NMR), we show that short ENM-BD simulations can reproduce the residue fluctuation profiles and the dominant protein motions with remarkable accuracy. Correlations of residue fluctuations often exceeded 80%, with many proteins showing agreements above 90%, while Principal Component Analysis (PCA) revealed strong correspondences between the essential motions in the ENM-BD simulations and those observed in MD and NMR ensembles. Following CG-to-all-atom reconstruction and short energy minimization, the ENM-BD conformers also achieve high stereochemical quality, which enables their usage in downstream atomistic applications. Finally, we show that the stochastic dynamics emerging from BD trajectories closely aligns with the normal modes encoded in the underlying elastic network, highlighting that most of the intrinsic dynamics of proteins is already embedded in their native 3D topology.
Dargis, R.; Arya, G.
Show abstract
DNA often functions in heterogeneous molecular systems containing proteins, lipids, polymers, and other materials. All-atom molecular dynamics simulations can be used to study DNA in these multicomponent systems, but computational cost limits the accessible system sizes and time scales. Coarse-grained models extend these scales, but existing DNA models are generally not designed for interactions with a broad range of other molecular species. To fill this gap, we develop a coarse-grained model of DNA designed for use with the Martini 3 force field. The model was parameterized through an iterative Bayesian optimization workflow, which used a scaled Wasserstein metric to compare distributions of local geometrical features and global structure from coarse-grained simulations against all-atom reference simulations. The optimized model captures key structural and mechanical properties of single- and double-stranded DNA across varying strand lengths and ionic conditions, while retaining compatibility with the broader Martini 3 ecosystem. This compatibility enables DNA to be integrated with a broad range of molecular systems, as we illustrate through simulations of double-stranded DNA bound to a transcription factor, cholesterol-tagged DNA duplex interacting with a lipid bilayer, a crossover-containing DNA nanostructure, and single-stranded DNA adsorbing onto graphene. Together, these results establish a transferable coarse-grained model of DNA for simulations of heterogeneous biomolecular and engineered systems.
Perez-Segura, C.; Scott, L. W.; Zlotnick, A.; Hadden-Perilla, J. A.
Show abstract
Protein function often depends on ligand binding pockets that fluctuate among conformational states, altering their size, shape, topology, and accessibility, yet quantitative comparison of these dynamic cavities remains challenging because their boundaries are often inherently ambiguous. The measure volinterior algorithm uses fuzzy-boundary detection to characterize enclosed molecular spaces; here, the hepatitis B virus (HBV) capsid assembly modulator (CAM) binding site is used as a model system to develop and validate a practical workflow for applying the method to dynamic protein binding pockets. The resulting methodology provides practical guidance for parameter selection and evaluation, establishes a standardized protocol for quantitative characterization of the HBV CAM pocket, and demonstrates robust, reproducible performance across conformational ensembles derived from molecular dynamics (MD) simulations. More broadly, this work provides a reproducible strategy for adapting measure volinterior to other dynamic binding pockets, enabling consistent comparison of pocket geometry among independent structural studies. Graphical Abstract O_FIG O_LINKSMALLFIG WIDTH=200 HEIGHT=105 SRC="FIGDIR/small/743403v1_ufig1.gif" ALT="Figure 1"> View larger version (30K): org.highwire.dtl.DTLVardef@d460b5org.highwire.dtl.DTLVardef@11915c2org.highwire.dtl.DTLVardef@1e37524org.highwire.dtl.DTLVardef@1fc1b7_HPS_FORMAT_FIGEXP M_FIG C_FIG
ADUPA, V.; Polet, J. D.; Dekker, M.; Onck, P. R.
Show abstract
Polyglutamine (polyQ) aggregation plays a central role in several neurodegenerative diseases, including Huntington's disease. DNAJB6b, a molecular chaperone involved in protein quality control, is known to efficiently suppress polyQ aggregation, but its anti-aggregation mechanism remains unclear. In this work we investigate the interaction between DNAJB6b and the polyQ region (Q48) of mutant Huntingtin Exon 1 (mHttEx1) using a custom-built coarse-grained molecular dynamics model. The model incorporates a two-bead-per-amino-acid representation with hydrogen bonding (termed 2BPA-HB), and is calibrated against all-atom molecular dynamics data in terms of geometry, hydrophobicity, and hydrogen bonding. The model reproduces the tertiary structure of DNAJB6b and its interactions with Q48, and reveals an inverse correlation between DNAJB6b concentration and Q48 aggregation propensity. Our simulations show that DNAJB6b co-condensates with polyQ molecules, thereby shielding the polyQ from forming the intermolecular hydrogen bonds necessary for amyloid formation. The 2BPA-HB CGMD model en- ables efficient exploration of DNAJB6b conformations, supporting future studies of chaperone-mediated aggregation suppression and therapeutic development.
Wu, D.; Wang, T.
Show abstract
The fidelity of molecular dynamics (MD) simulations fundamentally depends on the quality and coverage of the ab initio data used to parameterize the underlying force field, yet the role of side-chain conformational space remains insufficiently explored. In this study, we systematically investigate how comprehensive ab initio sampling of dipeptide conformations--specifically targeting side-chain degrees of freedom--impacts force field accuracy and MD simulation predictive power. We present the Protein Unit Dataset Plus (PUD+), a 40-million-conformation quantum mechanical dataset featuring unprecedented coverage of both backbone and side-chain conformational space. Machine learning force fields trained on PUD+ and integrated into AI2BMD simulations demonstrate superior energy and force prediction accuracy, capturing high-fidelity protein folding dynamics and the conformational flexibility of long-side-chain systems. Furthermore, leveraging PUD+ to reparameterize the CMAP term of the classical ff19SB force field markedly improves the description of intrinsically disordered protein (IDP) dynamics and IDP-ligand binding. Collectively, these results demonstrate that ab initio sampling of dipeptide side-chain conformations enables high-fidelity modeling of protein dynamics across both AI-driven and classical simulation paradigms.
Ardizzone, C.; Stegani, B.; Bacic Toplek, F.; Gianni, S.; Capelli, R.; Camilloni, C.
Show abstract
Accurately predicting how mutations alter protein-peptide binding remains challenging for molecular simulations because both conformational sampling and binding kinetics are computationally demanding. Here, we investigate whether multi-eGO, a hybrid transferable/structure-based atomistic-resolution model previously developed and validated for protein-small molecule interactions, can be transferred to protein-peptide binding without peptide-specific retraining. Using the PDZ2 domain of protein tyrosine phosphatase basophil-like in complex with the peptide EQVTAV as a benchmark, we first show that multi-eGO reproduces the structural dynamics of PDZ2 and the equilibrium binding thermodynamics of the wild-type complex. The model substantially accelerates both binding and unbinding relative to experiment but accurately preserves the resulting equilibrium dissociation constant. We then introduce conservative mutations in PDZ2 and in the peptide and evaluate their effects on binding without repeating the computationally expensive training procedure. Multi-eGO reproduces the experimentally observed changes in equilibrium dissociation constants, with strong agreement across PDZ2 mutants and moderate agreement when the peptide is also mutated. In contrast, the individual association and dissociation rate constants show substantially weaker agreement with experiment. The results indicate that the simplified energy landscape of multi-eGO limits the quantitative prediction of absolute kinetics while preserving thermodynamic information relevant to relative binding affinity. These findings establish multi-eGO as a computationally efficient approach for protein-peptide recognition and for predicting and rank-ordering the effects of conservative mutations on binding affinity.
Toplek, F. B.; Borges-Araujo, L.; Lindorff-Larsen, K.; Everaers, R.; Souza, P. C. T.; Morozova, T. I.
Show abstract
Biomolecular condensates formed by intrinsically disordered proteins require molecular models that accurately describe proteins in both dilute solution and condensed phases. Explicit-solvent coarse-grained models offer an attractive balance between chemical resolution and computational efficiency. Yet, it remains unclear whether improving dilute-state properties is sufficient to obtain an accurate description of condensates. Here, we address this question by introducing minimal modifications to the Martini 3 force field that combine recent advances in bonded interactions with refined protein-water interactions and strengthened glycine self-interactions, while preserving the underlying chemical transferability of the model. The resulting model substantially improves the description of single-chain conformations across a diverse benchmark of disordered proteins. We then investigate phase separation of the well-characterized low-complexity domain of heterogeneous nuclear ribonucleoprotein A1 and its sequence variants. The model reproduces several key physicochemical properties of biomolecular condensates, including chain expansion in the dense phase, sequence-dependent intermolecular contacts, protein diffusion and its relation to single-chain dimensions, and hydration, while revealing quantitative limitations in condensate density, phase equilibria, and ion partitioning. Our results show that improving dilute-state behaviour translates into a better description of condensed-phase properties, including condensate density, but is not sufficient to quantitatively reproduce the equilibrium between the dilute and dense phases.
Kudo, T.; Ekimoto, T.; Yamane, T.; Ikeguchi, M.
Show abstract
Many functional RNA motifs adopt structures that deviate from the canonical A-form helix and are emerging targets for RNA-directed therapeutics. The microtubule-associated protein tau (MAPT) A-bulge motif (5'-GCAGU/5'-ACGU) is one such motif. Because its structure is stabilized by a delicate balance of local interactions, its accurate modeling remains a major challenge for molecular dynamics (MD) simulations. The experimentally determined nuclear magnetic resonance (NMR) structure of the MAPT A-bulge motif provides a stringent test of whether RNA force fields can accurately reproduce the experimentally observed conformation. Most current AMBER-family RNA force-field models have incorrectly favored a non-native base-triple state of the MAPT A-bulge motif over the experimentally observed stacked state. Structural comparison of the stacked and base-triple conformations revealed that overly favorable NH-N hydrogen bonds between the bulged adenosine and an adjacent Watson-Crick base pair were the primary source of this imbalance. We developed gHBfix-18Ab, an 18-component hydrogen-bond correction that distinguishes NH and NH2; donors. gHBfix-18Ab was combined with the previously developed OL3CP and NBfix0BPh corrections to generate the composite model gHBfix-18Ab*. This model restored the experimentally observed stacked state as the global minimum in the calculated free-energy profile and improved agreement with NMR-derived distance data for the A-bulge region. Importantly, gHBfix-18Ab* did not produce marked structural destabilization of the cUUCGg tetraloop, a widely used benchmark for RNA force-field validation, suggesting that the refinement preserves the stability of the unrelated RNA motif. These results demonstrate that targeted refinement of hydrogen-bond interactions provides a practical strategy for systematic improvement of RNA force fields toward more accurate modeling of noncanonical RNA motifs.
Yasuda, I.; Tesei, G.; Yamamoto, E.; Yasuoka, K.; Lindorff-Larsen, K.
Show abstract
Biomolecular condensates function as membraneless compartments, and some protein condensates can selectively concentrate single-stranded nucleic acids while excluding double-stranded nucleic acids. Understanding how nucleic acid structure affects partitioning into condensates has important implications for nucleic acid activity and function within condensates. Here, we present a set of coarse-grained two-bead-per-nucleotide models for simulations of double-stranded RNA and DNA in the CALVADOS framework. Our models separately represent the backbone and base, and maintain the helical structures using an elastic network potential tuned to capture chain stiffness. For dsRNA, the base stickiness was tuned using experimental data on differential partitioning of single- and double-stranded RNA into Ddx4N1 condensates in order to account for reduced base accessibility upon duplex formation. This RNA structural selectivity varied with the balance of electrostatic and non-electrostatic interactions, as revealed by simulations of condensates of the CAPRIN1 disordered region at varying ionic concentrations and with an R-to-K sequence variant. Finally, we developed parameters for double-stranded DNA using a similar approach. We envision that the CALVADOS models for double-stranded RNA and DNA will be useful for studying co-condensates of proteins and structured nucleic acids.
Yang, Y.
Show abstract
Gaussian accelerated molecular dynamics (GaMD) enhances conformational sampling by adding a smooth boost potential without requiring predefined collective variables, but an engine-integrated implementation has not been available in GROMACS. Here, we implement total-, dihedral-, and dual-boost GaMD in GROMACS 2025.4, including staged energy-statistics collection, GPU-based bias evaluation and force scaling, restart support, and outputs required for cumulant-based free-energy reweighting. The implementation was evaluated using four benchmark systems spanning conformational free energies, protein folding, and ligand recognition. For alanine dipeptide, a reweighted 100 ns GaMD trajectory recovered the major free-energy basins and rotational barriers in overall agreement with a 1000 ns conventional MD simulation. For chignolin and TC5b, all three independent trajectories for each system sampled native-like folded states from extended conformations within 300 ns and 1 s, respectively; the best TC5b structure had a minimum backbone RMSD of 0.03 nm from the experimental structure. In the benzene-T4 lysozyme system, two of five independent 500 ns trajectories captured both ligand binding and dissociation, yielding a bound pose with a minimum ligand RMSD of 0.06 nm from the crystal structure. Across all four systems, the boost-potential distributions were approximately Gaussian, and second-order cumulant reweighting resolved the expected conformational and binding free-energy basins. These results demonstrate that GROMACS-GaMD provides a practical, GPU-enabled, collective-variable-free enhanced-sampling framework for biomolecular free-energy calculations, protein folding, and ligand-binding studies.
Ojha, A. A.; Huber, G.; Dutta, S.; Hanson, S. M.
Show abstract
Drug-target association and dissociation rates often determine in vivo efficacy more than affinity alone. The association rate constant, however, is computationally challenging to predict since the productive encounter is a rare event in a large translational and orientational search space, which lies beyond the reach of conventional atomistic simulations. We present PySTARC (Python Simulation Toolkit for Association Rate Constants), a GPU-accelerated Brownian dynamics engine for estimating bimolecular association rate constants. PySTARC is a Python reimplementation of the BrownDye engine that converges even small reaction probabilities on a single GPU, at a throughput that would otherwise require a large CPU cluster. The current framework resolves reactions between integration steps with a closed-form Brownian bridge, models the internal flexibility of the solute through a coarse-grained bead chain, distributes trajectories across multiple GPUs, checkpoints long runs, monitors convergence, and automates the entire workflow from input structures to rate estimates. PySTARC is validated against protein-ligand and protein-protein complexes spanning five orders of magnitude in the association rate constant, with most estimates within one order of magnitude of experiment.
Nadeem, H.; Kleiman, D. E.; Shukla, D.
Show abstract
Adaptive sampling accelerates the exploration of conformational space in molecular dynamics (MD) simulations by repeatedly analyzing the accumulated trajectories and seeding a new round of simulations from informative configurations. A growing collection of adaptive sampling policies has been proposed, each built around a particular notion of what makes a configuration informative, yet these methods are scattered across separate and often incompatible implementations, which complicates their systematic comparison and their combined use in meta adaptive sampling schemes. Here, we present AdaptivePy, a compact and extensible Python framework that implements nine seed-selection policies behind a single configuration-driven interface, spanning simple population-based baselines, several established machine-learning and geometry-based methods, and two ensemble or meta sampling policies introduced in this work. We show that the shared implementation reproduces the characteristic selection behavior of each policy on a series of analytic benchmark landscapes. We also introduce a new adaptive sampling scheme that employs TS-DAR, a deep learning framework originally designed to identify transition states, into an acquisition criterion that drives the discovery of an entire multi-basin landscape starting from a single basin. We further demonstrate that the common interface enables meta adaptive sampling policies, which aggregate the rankings of several policies into a single set of seeds. AdaptivePy thereby provides a unified testbed for the adoption, benchmarking, and continued development of adaptive sampling methods for biomolecular MD simulations.
Noske, J.; Janzen, M.; Lepoivre, T.; Höcker, B.
Show abstract
Accurate prediction of mutation-induced changes in protein-protein affinity remains a central challenge in computational biophysics and protein engineering. Here we present a physics-based scoring method for predicting the effects of single amino acid substitutions on protein-protein and protein-peptide binding affinity. The method is based on MMPBSA energy calculations and implemented in a user-friendly, publicly available workflow (https://github.com/Hoecker-Lab/mmpbsa_scoring). To derive the scoring function, MMPBSA energy terms were fitted using a large and diverse dataset of experimentally measured binding affinity changes upon point mutation. The terms were then evaluated on independent systems. Specifically, we assessed the ability to predict mutational effects at the interface of SARS-CoV-2 spike protein with its cellular receptor, as well as computationally predicted Armadillo-repeat protein-peptide complexes. The approach performs well across diverse targets and shows promising predictive performance on both experimentally determined and computationally predicted structures. It matches or exceeds the accuracy of established physics-based protein design software and state-of-the-art deep learning-based scoring functions when applied to data outside their training domain. The method provides a computationally efficient and interpretable framework for prioritizing interface mutations in protein engineering and binding affinity optimization.
Prajapati, S.; Mathur, N.; Chakraborty, D.; Srivastava, A.
Show abstract
We present SOP-Multi-{pi}, an enhanced version of the two-bead-per-amino-acid resolution Self-Organized Polymer (SOP) model designed to investigate the conformational dynamics of Intrinsically Disordered Proteins (IDPs) and their phase behavior in biomolecular condensates. The original SOP-IDP and SOP-Multi models work well for chains in good solvents that predominantly adopt expanded conformations. However, they fail to reproduce the single-chain properties and phase-separation behavior of IDPs such as hnRNPA1 and its variants, which have a distinct pattern of aromatic residues. We therefore reparameterized the SOP-IDP force field to explicitly capture {pi}-{pi} and cation-{pi} interactions at the coarse-grained (CG) resolution. These interactions are weak but important for biomolecular recognition and for driving and stabilizing biomolecular condensates. Dimensionless rescaling prefactors ({lambda}) were applied exclusively to the aromatic {pi}-{pi} (Phe, Tyr, Trp) and cation-{pi} (Arg/Lys with Phe, Tyr, Trp) pairs. They were optimized iteratively against experimental radii of gyration of the hnRNPA1 low-complexity domain (A1-LCD) and its 16 aromatic and charge mutants. This restored a physically meaningful Tyr > Phe interaction hierarchy, and the Pearson correlation for Rg improved from r = 0.07 (SOP-IDP) to r = 0.83 (RMSD = 0.08 nm). The same parameters, without further tuning, gave r = 0.93 for 34 additional IDPs (24-198 residues) and reproduced Flory scaling ( {nu} ~ 0.588). We back-mapped microsecond-scale CG ensembles to all-atom representations and computed SAXS and Kratky profiles for direct comparison with experiments on IDPs and on the multidomain proteins hnRNPA1, TDP-43, G3BP1, and FUS. In direct-coexistence slab simulations, the model captures the phase behavior of A1-LCD condensates and their variants, with critical temperatures that correlate strongly with experiment (r = 0.88) but are somewhat underestimated in absolute terms. The model is implemented in LAMMPS, and a detailed user manual is provided for setting up single-chain and multi-chain condensate simulations of full-length multidomain proteins.
Li, D.; Ken, M.
Show abstract
Accurate determination of RNA conformational ensembles is essential for understanding RNA function and advancing RNA-targeted drug discovery, yet lowly-populated alternative states remain difficult to resolve with atomistic detail. A central constraint is that experimental refinement can only select conformations already present in the starting library, making library generation the limiting step. Using the HIV-1 trans-activation response element (TAR) as a model system, we benchmarked conventional MD (cMD) against the enhanced sampling methods Gaussian-accelerated MD (GaMD), replica-exchange Gaussian-accelerated MD (Rex-GaMD), replica-exchange with solute tempering (REST2), and temperature replica-exchange MD (T-REMD), as well as the structure-prediction based methods FARFAR2 and AlphaFold 3. Each library was refined against experimental residual dipolar couplings (RDC) and validated independently using ensemble-averaged QM/MM chemical shifts. We showed that T-REMD produced the most accurate ensemble by both measures, and its advantage tracked with broader, more continuous coverage of the interhelical conformational landscape. Broad temperature-range T-REMD also sampled conformations resembling excited state 1 (ES1) and the U23-A27-U38 base-triple, without requiring these states to be specified during library generation. More accurate ensembles further improved coverage of experimentally observed ligand-bound TAR conformations and enhanced ensemble-based virtual screening, linking structural accuracy to functional utility. The same workflow applied to the preQ1 class I riboswitch and the UUCG tetraloop improved agreement with experimental data in both cases. Together, these results establish replica-exchange enhanced sampling, particularly T-REMD, as an effective strategy for constructing experimentally validated RNA ensembles and accessing conformations corresponding to rare functional substates.
Bhandari, T. R.; Kremer, K.; Girard, M.
Show abstract
Simulations have emerged as a pillar in biophysics to understand behavior at the molecular scale, in particular for proteins such as intrinsically disordered proteins. These often show transient folding with long-lived states that are challenging to efficiently sample using conventional simulations. Here, we show that Hamiltonian replica-exchange at the level of hydrogen bonds can be used to systematically explore the conformational landscape of proteins. Investigating the IM30 protein, a partially alpha helical protein, under the PLUM coarse-grain force-field, we show that this approach straightforwardly allows characterization of local helix stability, which would be otherwise hard to obtain. For instance, we show that some regions of the conserved coiled-coil structure are less stable than others. By combining the different replicas with principal component analysis, we are able to dissect the conformational landscape of the IM30 protein. Taken together, this indicates that enhanced sampling remains an important part of simulations, allowing insights that cannot be easily obtained in traditional simulations.
Li, Y.; Lu, K.
Show abstract
Conformation ensembles of biomolecules provide the basis for understanding structural transformations and drug design. Deep-learning generative models have advanced protein and small molecule ensemble generation, while RNA-Protein complexes remain unaddressed due to the chemical heterogeneity, limited dataset size and the different flexibility scales of RNA and protein components. We present RPDynaFlow, a flow-matching model to generate conformation ensembles of RNA-protein complexes, trained on 600 ns trajectories of molecular dynamics(MD) simulation. The results show our model extends the sampling range of the phase space compared to MD simulation, which couldbe treated as a rapid and efficient complement to MD trajectoriesfor studying RNA-protein interactions.